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ABSTRACT 

Context. Deuterated ions, especially H 2 D + and N2D + , are abundant in cold (~ 10 K), dense (~ 10 5 cm -3 ) regions, in which CO is frozen out 



O ' onto dust grains. In such environments, the N 2 D + /N 2 H + ratio can exceed the elemental abundance ratio of D/H by a factor of = 10 
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Aims. We use the deuterium fractionation to investigate the evolutionary state of Class protostars. In particular, we expect the N 2 D + /N 2 H + 



^2 ' ratio to decrease as temperature (a sign of the evolution of the protostar) increases. 

, Methods. We observed N 2 H + 1-0, N 2 D + 1-0, 2-1 and 3-2, C ls O 1-0 and HCO + 3-2 in a sample of 20 Class and borderline Class 0/1 protostars. 
We determined the deuteration fraction and searched for correlations between the N 2 D + /N 2 H + ratio and well-established evolutionary tracers, 
i-H ■ such as T Dust and the CO depletion factor. In addition, we compared the observational result with a chemical model. 

Results. In our protostellar sample, the N 2 H + 1-0 optical depths are significantly lower than those found in prestellar cores, but the NtH + 
column densities are comparable, which can be explained by the higher temperature and larger line width in protostellar cores. The deuterium 
fractionation of N 2 H + in protostellar cores is also similar to that in prestellar cores. We found a clear correlation between the N 2 D + /N 2 H + ratio 
and evolutionary tracers. As expected, the coolest, i.e. the youngest, objects show the largest deuterium fractionation. Furthermore, we find 
that sources with a high N 2 D + /N 2 H + ratio show clear indication for infall (e.g. Sv < 0). With decreasing deuterium fraction the infall signature 
ON disappears and dv tends to be positive for the most evolved objects. The deuterium fractionation of other molecules deviates clearly from that 
CD of N 2 H + . The DCO + /HCO + ratio stays low at all evolutionary stages, whereas the NH 2 D/NH 3 ratio is > 0. 15 even in the most evolved objects. 
Conclusions. The N 2 D + /N 2 H + ratio is known to trace the evolution of prestellar cores; we show that this ratio can be used to trace core 
evolution even after star formation. Protostars with an N 2 D + /N 2 H + ratio above 0.15 are in a stage shortly after the beginning of collapse. Later 
on, deuterium fractionation decreases until it reaches a value of ~ 0.03 at the Class 0/1 borderline. 
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1. Introduction tween the bolometric luminosity and the distance normalized 

1.3 mm flux (Saraceno et al. |19961 l. All these evolutionary in- 

Early attempts to describe an evolutionary sequence of proto- m bafjed Qn ^ change of ^ ^ distdbu _ 

stars classified them according to their near- and mid-infrared . • /( , cr ,^ f . u . . A .. , . • 

2, tion (SED) of the protostar. As it heats up, an increasing frac- 

spectra. Lada & Wilking <T5P» divided a sample of protostel- ^ Qf ^ radiation ^ emined ^ shorter wavel ^ leading 

lar sources in the p Ophiuchi cloud into three classes (Class I to tQ an increase of Lbql and Tbol whh ^ 

III with progressive evolutionary stage). Andre et al. ( 1993) in- _ , , . , , , . , . 

, . „ , .... . . Several theoretical models connect physical properties of 

troduced the Class stage. The strong submillimetre emission , • , , • • • • , „ , 

„ _ „ , . , . , • , , • , ■ • the protostars with the time since gravitational collapse (e.g. 

of Class objects relative to their bolometric luminosity sug- _ ., i ,„„ i . „„ ol TT . . . 

, j. , ■ „ i • i , Smith 1998 Myers et al. 1998). However, although the evo- 



gests that the mass of the circumstellar envelope is larger than 
the mass of the central star. Class objects are the youngest 
among the protostars. In the subsequent years, more sensi- 



lutionary stages determined by different models is the same 
(i.e. the time sequence is the same), the absolute ages vary sig- 



, , , ,. , mficantly. The absolute age for Class 0/1 borderline objects, 

tive indicators for protostar evolution have been established, . , . , ni . . , 

, ... „„ „ . ..rrrnrn , tor example, varies between 10 and a few times 10 years 

such as the bolometric temperture (Myers & Ladd 1993), the , . f i .„„^ t T , . , ... ... 

T „ . ,. . , , . . . . . (Froebnch 2005). In this paper, we show that it is possible 

L-BOL/L-smm ratio (Andre et al. 1 1993b and the correlation be- , . „ . . „. 

to use the deuterium fractionation of N2H as an evolutionary 
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freeze out onto dust grains (e.g. Willacy et al. [T998 , Caselli et 
al. [19991 Kramer et al. [19^91 Bergin et al. [20021 Redman et 
al. 120021 Tafalla et al. 121)021 120041 . By comparison, molecules 
such as N2H + and NH3 suffer less from depletion. The persis- 
tence of nitrogen-bearing species is likely due to the fact that 
N2, the progenitor of these species, remains in the gas phase 
longer than CO, despite the fact that recent laboratory mea- 
surements have found that the binding energy, as well as the 
sticking coefficient of N2 and CO on dust, are quite similar 
(Oberg et al. 2005, Bisschop et al. I20061 >. Recent theoretical 
work suggests that the resistance of nitrogen-bearing species to 
depletion may be due to the fact, that a significant fraction of 
interstellar nitrogen is in atomic (N) rather than molecular (N2) 
form (Maret et al. 2006, Flower et al. 2006). Atomic nitrogen 
has a smaller binding energy than molecular nitrogen and CO. 

One consequence of the low temperatures and freeze-out 
of CO in pre- and protostellar cores is that the abundance 
of deuterated molecules is greatly enhanced. In such regions, 
column density ratios for N 2 D + /N 2 H + of 0.24 (Caselli et 
al. I2002bl l. D 2 CO/H 2 CO 0.01-0.1 (Bacmann et al. 120051 . 
NH2D/NH3 up to 0.33 (Hatchell 120031 and DCO + /HCO + in 
the order of 0.01 (J0rgensen et al. 2004) are found. These ratios 
are three to four orders of magnitudes larger than the average 
ISM D/H ratio of 1.510 5 (Oliveira et al. 120031 . Models show 
that the abundance of H2D + , a fundamental molecule in deu- 
terium chemistry, rises steeply as temperature decreases (e.g. 
Flower et al. 120041 Roberts & Millar [2000). Consequently, the 
deuterium fractionation of many other molecules also increases 
(Millar etaLEpOOJ. 

In a recent study, Crapsi et al. (2005) investigated the 
N2D + /N2H + ratio in 31 low-mass starless cores and found 
ratios in the range of a few percent to 0.44. They also found 
a tight correlation between the N2D + /N2H + ratio and the CO 
depletion factor; the ratio is greater where CO is depleted. 
Furthermore, they attempted to find correlations between the 
N2D + /N2H + ratio and various chemical and dynamical evolu- 
tionary indicators (e.g. core density, molecular column density, 
line width, and line asymmetry). Although the dependencies 
on particular indicators are not always clear, they identified a 
recognizable trend. The N2D + /N2H + ratio was generally found 
to increase as a cloud core evolves towards the moment of pro- 
tostellar collapse. Another results from their work is that the 
trends seen in the subsample of Taurus cores were more homo- 
geneous than trends in the total sample. This homogeneity led 
them to speculate that a core's external environment plays an 
important role in its evolution. In this work we follow the trend 
of the N2D + /N2H + ratio to the subsequent, early protostellar 
stages. 

One recent study that also sought to determine the behavior 
of deuterated molecules in protostellar environments is from 
Roberts & Millar J2007l . They determined the N 2 D + /N 2 H + ra- 
tio in five low mass protostellar cores, using the University 
of Arizona 12m telescope (HPBW of 70" for N 2 H + 1- 
0). Furthermore they measured the D2CO/H2CO and the 
HDCO/H2CO ratio in these cores. The HDCO/H 2 CO ratio was 
similar in most cores (~ 0.06), with the exception of the ob- 
ject with the highest N 2 D + /N 2 H + ratio, where HDCO/H 2 CO 
was lowest. For four sources, the N2D + /N2H + ratio appeared 



anti-correlated with the D 2 CO/H 2 CO ratio. Only in L 1527 
were both N2D + and D2CO abundances anomalously low. They 
did not find a correlation between deuterium fractionation and 
bolometric temperature for any of the observed molecules. 
However, due to the small number of protostars, this work can- 
not rule out any trend of the deuterium fractionation in proto- 
stellar cores. 

We observed 20 Class protostars in N2H + and N2D + , us- 
ing the IRAM 30m telescope (HPBW of 27" for N 2 H + 1-0, 
see section©. The goal of our study is to investigate the trend 
of deuterium fractionation after the moment of star formation. 
Spectra and line parameters are presented in Sect. [5] We expect 
that the N2D + /N2H + ratio (determined in section |4~T1 ) should 
decrease during protostellar evolution due to the internal heat- 
ing of the core. As the core environment warms up, H2D + 
and its more highly deuterated isotopologues are destroyed. 
Consequently, all other deuterated molecules formed by re- 
actions with H2D + (such as N2D + ) are destroyed as well. In 
Sect. l4~2l |44l the N 2 D + /N 2 H + ratio is compared with evolu- 
tionary indicators, namely dust temperature, CO depletion fac- 
tor, and Lbol/Fi.3- Possible correlations with kinematical pa- 
rameters of the gas are investigated in section 1431 In addition, 
we compare the N2D + /N2H + ratio with the deuterium fraction- 
ation of other molecules ( section l4~6l . Finally, we compare our 
results with the predictions of a simple chemical model (sec- 
tion [5]). 



2. Observations 

2.1. Source sample 

Our survey includes 20 Class protostellar cores. Tab. Q] lists 
the objects, their coordinates and their estimated distances. The 
sources were selected in terms of being a.) well known (i.e. 
many complementary data can be found in the literature) and 
b.) representing a range of stages of Class sources, from 
very young objects (e.g. HH 211, which is one of the youngest 
protostars according to Froebrich 2005 ) to borderline objects 
(Class 0/1), such as L1455 Al (Froebrich 120051 and Barnard 5 
IRS 1 (Velusamy & Langer fT9"9"8l . For one object, L 1527, the 
classification is ambiguous. The continuum data denote it as a 
Class object, but an embedded source is detected at a wave- 
length < 5 ymx (Froebrich 120051 ). indicating that this protostar 
is already more evolved. Therefore, this source is also classi- 
fied as a Class 0/1 borderline object. Most of the objects (13) 
are located in the Perseus molecular cloud, and all are at a dis- 
tance closer than 500 pc. For the sources located in Perseus, 
we assumed a distances of 220 pc (Cernis 1990). Other authors 
(e.g. Herbig & Jones 1983 ) give distances larger than 300 pc. 
Eleven objects belong to a sample of sources catalogued by 
Froebrich (2005 ), and for 17 sources, kinetic gas temperatures 
are known from NH3 observations (Jijina et al. 1999). However, 
the NH3 observations vary in spatial resolution between 40" 
and 80". To avoid systematic errors, introduced by the differ- 
ent beam sizes, we use these kinetic temperatures here only as 
an estimate of the gas temperature. 
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Table 1. Class objects observed for the survey. 



Source 


D A 

KA 




Distance 




(MUUU.U) 


(JZUUU.U) 


LpcJ 


L 1448 IRS 2 


03:25:22.4 


30:45:12 


220 1 


L 1448 IRS 3 


03:25:36.4 


30:45:20 


220 1 


L 1448 C 


03:25:38.8 


30:44:05 


220 1 


L 1455 Al 


03:27:42.1 


30:12:43 


220 1 


L 1455 A2 


03:27:49.8 


30:11:42 


220 1 


Per4S 


03:29:17.9 


31:27:30 


220' 


Per4N 


03:29:22.7 


31:33:30 


220 1 


Per 5 


03:29:51.6 


31:39:03 


220 1 


Per 6 


03:30:14.8 


30:23:48 


220' 


IRAS 03282+3035 


03:31:21.0 


30:45:30 


220 1 


Barnard 1-b 


03:33:21.3 


31:07:35 


220 1 


HH211 


03:43:56.8 


32:00:50 


220 1 


Barnard 5 IRS 1 


03:47:41.6 


32:51:42 


220 1 


L 1527 


04:39:53.5 


26:03:05 


140 2 


L483 


18:17:29.8 


-04:39:38 


200 2 


SMM5 


18:29:51.2 


01:16:40 


310 3 


L723 


19:17:53.1 


19:12:16 


300 2 


L673 A 


19:20:25.8 


11:19:51 


300 2 


B 335 


19:37:01.4 


07:34:06 


250 2 


L 1157 


20:39:06.5 


68:02:13 


440 4 



'Cernis fT990t . 2 Jijina et al. {[999 1, 3 Froebrich J20051 . 4 Gueth et 
al. ( [20031 . 

2.2. Molecular line observations 

The observations were carried out in January 2006 using the 
IRAM-30m telescope. In Tab. [2] the molecular line transi- 
tions, their frequencies, and their corresponding half power 
beam widths (HPBW) are given. Two different receiver set-ups 
were used. The first set-up used frequency switching to observe 
C I8 1-0, N 2 D + 2-1, N 2 D + 3-2 and HCO + 3-2 simultaneously. 
The offsets were 10 (C lg O 1-0), 15 (N 2 D + 2-1), 25 (N 2 D + 3-2) 
and 28 (HCO + 3-2) MHz, respectively. The second set-up used 
an external local oscillator to observe the N 2 D + 1-0 line. In 
this mode, frequency switching was not available, so position 
switching was used instead. The N 2 D + 1-0 line was observed 
simultaneously with the N 2 H + 1-0 and N 2 D + 3-2 lines. As a 
back end, the VESPA autocorrelator was used, which was ad- 
justed to a resolution of about 0.05 km/s at each frequency. 

Table 2. Observed emission lines. 



Line Frequency [GHz] HPBW ["] 



N 2 D + 1-0 77.109 32 

N 2 H + 1-0 93.173 27 

C 18 1-0 109.782 22 

N 2 D + 2-1 154.217 16 

N 2 D + 3-2 231.322 11 

HCO+ 3-2 267.558 10 



During observations, the atmospheric r(225 GHz) was 
about 0.25, which is a typical opacity for winter weather con- 



ditions at the IRAM site. The corresponding system tempera- 
tures were between 150 K (N 2 H + 1-0) and 800 K (HCO + 3- 
2). Pointing was checked every two hours, and the pointing 
accuracy was better than 3". Antenna temperatures and main 
beam temperatures are related by the main beam efficiency, 
?7mb, which is between 0.8 and 0.42, depending on the observ- 
ing frequency^. 

3. Results 

3.1. N 2 H + &N 2 D+ 

The results of the N 2 H + and N 2 D + observations are summa- 
rized in Tab. [3] Fig. [TJ shows the three N 2 D + transitions ob- 
served toward HH 211, one of the few sources (total of 6, see 
TableO where all lines are detected. The N 2 H + 1-0, N 2 D + 2-1, 
C ls O 1-0 and HCO + 3-2 spectra of all 20 sources are displayed 
in Fig. [2] 
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- N 2 D + 3-2 




H — i — i — i — — i — i — i — i — — i — ^ 

: N 2 D + 2-1 


— i — i — — i — i — i — i — — i — i — i — ^ 

— 1 — 1 — — 1 — 1 — 1 — 1 — — 1 — 1 — 1 — M 


~—\ — I — I — I — — I — I — I — I — — I — ^ 

N 2 D + 1-0 





-5 5 10 15 20 

Velocity (km/s) 

Fig.l. Spectra of N 2 D + 1-0, N 2 D + 2-1 and N 2 D + 3-2 toward 
HH 211, one of the youngest sources in our sample. 

The N 2 H + 1-0 line was detected in all sources of our sam- 
ple. The total optical depth was determined by fitting all hy- 
perfine components of the line using the "METHOD HFS" of 

1 http://www.iram.fr/ 
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N 2 H (1-0) 




BARNARD1 — B 



A 



-10 1 

Velocity-V LSR [km/s] 



N 2 D + (2-1) 



i^p'wIj'Kj^w.-ASV.in 



BARNARD1 — B 



A>/^.^.«'-'<> i '..v^/-Y»''iWii 



^^yp^^^rfV jV^ 1 



j^Y^'»*l|V v ''|V^ / 'j'V 



-10 10 

Velocity-V LSR [km/s] 



c 18 o (1-0) 







-10 

Velocity-V LSR [km/s] 



HC0 + (3-2) 



^wVjW^jpyy^ ^ 



BARNARD1 — B 



A 



[_ PER6 



'''' ■• ^ 



Velocity-V LSR [km/s] 



Fig. 2. Spectra of N 2 H + 1-0, N 2 D + 2-1, C ls O 1-0 and HCO + 3-2 (from left to right) of all 20 protostars of our sample. To be able 
to plot all spectra in one plot we divided the intensities of HCO + by three. 



the CLASS progra The relative intensities of the hyperfine was estimated using the equation 
components given by Womack et al. ( 1992) were assumed. _ T 
With the exception of Per 6, the line optical depths could be Tmb = (J(Tex) ~ J ( Tbg ^ ' ^ ~ e ) 

determined with at least 3.5 cr in all objects (see Tab. |3}. T ex where T ex is the excitation temperature, T BG the cosmic back- 
ground temperature (2.7 K) and J(T) = [hv/k]/[exp(hv/kT) - 

1]. In most cases, the excitation temperature of N 2 H + is lower 

2 http://www.iram.fr/IRAMFR/GILDAS/ than the kinetic gas temperature listed in Jijina et al. ( 1999J) and 
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Table 3. Line parameters of the N2H + and N2D + lines. 







N2H + 


(1-0) 










N 2 D + 








Source 


w 


Av 


T 


T 

1 ex 


T 

iKin 


YI/ 


VVy=2-l 




T 


A ,7 

AV 2 -1 


T 

1 ex 




LJMcm/sJ 


[Km/Sj 






1 1. 1 


LKJcm/sJ 


[KKm/sJ 


[JSJOTI/SJ 


T O 1 


[Km/sJ 




L 1448 IRS 2 


12.93 


0.50 


11.6 


6.8 


- 


1.7 


1.10 


0.73 


< 0.1 


0.53 


- 


L 1448 IRS 3 


24.65 


0.93 


4.1 


11.4 


12.0 


4.0 


5.05 


3.22 


2.6 


0.59 


7.9 


L 1448 C 


12.48 


0.76 


5.3 


7.3 


13.7* 


< 1 


2.32 


1.38 


< 0.1 


0.72 


- 


L 1455 Al 


7.20 


0.74 


1.4 


9.2 


15.0 


< 1 


0.68 


<0.74 


< 0.1 


0.42 


- 


L 1455 A2 


8.60 


0.54 


3.0 


9.1 


11.0 


< 1 


1.03 


<0.73 


< 0.1 


0.33 


- 


Per 4 S 


5.08 


0.29 


3.5 


8.5 


9.0 


1.3 


1.43 


0.98 


< 0.1 


0.32 


- 


Per4N 


4.62 


0.31 


3.6 


7.6 


10.0 


1.3 


0.75 


<0.68 


< 0.1 


0.31 


- 


Per 5 


7.53 


0.38 


3.8 


8.7 


11.0 


< 1 


1.07 


1.08 


< 0.1 


0.36 


- 


Per 6 


4.64 


0.36 


<0.1 




11.0 


< 1 


< 0.18 


<0.69 








IRAS 03282 


9.53 


0.38 


5.3 


8.6 


11.3* 


< 1 


3.58 


1.47 


3.9 


0.34 


7.9 


Barnard 1-B 


19.05 


0.82 


8.6 


7.6 


12.0 


5.9 


7.48 


4.82 


4.9 


0.78 


6.7 


HH211 


9.18 


0.41 


6.2 


7.6 


12.8* 


2.3 


5.24 


2.08 


4.9 


0.39 


7.5 


Barnard 5 IRS 1 


8.59 


0.43 


6.4 


7.2 


10.8* 


< 1 


0.50 


0.65 


<0.1 


0.50 




L 1527 


3.91 


0.32 


5.6 


5.7 


10.8* 


< 1 


< 0.18 


<0.59 








L483 


14.05 


0.44 


12.1 


7.3 


10.0 


2.1 


1.75 


1.13 


2.8 


0.37 


5.5 


SMM5 


15.39 


0.70 


4.9 


8.9 


12.0 


< 1 


2.46 


1.43 


<0.1 


0.58 




L723 


5.70 


0.71 


2.4 


6.6 


11.0 


< 1 


0.34 


1.29 


<0.1 


0.32 




L 673 A 


7.34 


0.54 


2.9 


8.6 


<12.0 


< 1 


0.45 


0.63 


<0.1 


0.34 




B 335 


6.05 


0.40 


9.0 


5.7 


9.9 


< 1 


0.42 


0.42 


<0.1 


0.46 




L 1157 


8.38 


0.65 


2.1 


8.9 




< 1 


0.53 


1.00 


<0.1 


0.56 





The kinetic gas temperatures marked with an * are taken from Hatchell ( 2003 ). All other Tki n values are listed in Jijina et al. ( 1999 1. The errors 
on the integrated intensities of N 2 H + 1-0 and N 2 D + 2-1 are ~ 10% (calibration uncertainties), whereas the errors on the N 2 D + 1-0 and N 2 D + 3-2 
lines are dominated by the noise of the spectra (0.5 and 0.2 Kkm/s for J= 1-0 and J=3-2, respectively), t is the total optical depth, defined as the 
sum of the optical depths of the individual hyperfine components. Typical errors of r are ±0.5 and +0.2 for N 2 H + and N 2 D + 2-1 (where it is 
detected), respectively. The resulting error of T cx is ~ 0.2 K for both species. 



Hatchell (2003 ), most likely due to line of sight average volume 
densities of the emitting gas lower than the critical density of 
the observed transitions (^ 10 5 cm -3 for N2H + (l-0)). The fac- 
tor by which T ex is lower than T^, varies from source to source 
from close to 1 (L 1448 IRS 3) to 0.43 (HH 211). However, 
these excitation temperatures are typically 2-3 K higher than 
those in the prestellar cores observed by Crapsi et al. (J2005 ). 
Another contrast with the prestellar cores, is that none of the 
lines in our sample show highly non-Gaussian profiles. 

The N2H + 1-0 line widths found in our protostellar sample 
are significantly larger than the ones in prestellar cores. The av- 
erage line width for our sample is 0.61 kms" 1 , whereas Crapsi 
et al. d2003T > found an average line width of 0.26 kms 1 . The 
line widths reported by Roberts & Millar (2007 ) are slightly 
broader than ours (0.85 kms -1 ), but this difference might be 
due to the lower spatial resolution of their observations. This 
suggests that nonthermal motions are not only concentrated 
nearby the protostar, but that they pervade the whole associated 
envelope,increasing outward as commonly found in molecular 
cloud cores (e.g. Fuller & Myers 1992). 

N2D + was observed in three different transitions (J=l-0, 
J=2-l & J=3-2). All three lines were detected in six objects. 
The ground-state transition, N2D + 1-0, was detected in only 
seven sources, but the S/N ratio of four of these is too low 
to obtain a reliable r from hyperfine component fitting. Only 
in L 1448 IRS 3 and Barnard 1-b is the line strong enough; 
the HFS method yielded t = 1.3 + 0.4 and 0.9 + 0.2, re- 



spectively. The reason for the rarity of N2D + 1-0 detections 
is that the Einstein-A coefficient, i.e. the transition probability, 
increases with ~ J 3 , and thus the A-coefficient for the J=l-0 
line is ten times lower than the A-coefficient of the J=2-l line. 
Furthermore, the relatively low energies (E = J(J+ 1) • 1.85 K) 
and the higher degeneracies of higher rotational states lead to 
an efficient population of those states under conditions present 
in protostellar cores (n^ 10 5 cm -3 , T=10-20 K). 

The 2-1 line is the most frequently detected N2D + line. 
Only two objects, Per 6 and L 1527, do not show emission in 
J=2-l. These two objects also show no emission in the other 
N2D + lines. In most cases, the J=2-l lines are optically thin, so 
that t can be determined by fitting the hyperfine components in 
only five cases. 

The N2D + 3-2 line is detected in 15 out of the 20 cores. 
Because the satellite lines get weaker with increasing J, and 
due to the higher noise-level of most of these data, it was not 
possible to derive reliable optical depths. 

The integrated intensities of the N2H + 1-0 lines from the 
cores in our sample are comparable to the intensities found in 
prestellar cores (Crapsi et al. 2005 ), but the optical depths in the 
protostellar cores are typically lower by a factor of three. The 
integrated intensity and optical depth of N2D + 2-1 are approx- 
imately the same in pre- and protostellar cores, but N2D + 3- 
2 is significantly stronger in the latter (~ 0.3 Kkms -1 and 
~ 1.42 Kkms -1 , respectively). 
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3.2. C 18 & HCO + 

The C I8 1-0 and HCO + 3-2 lines were detected in all 20 
sources of our sample. The line parameters are listed in Tab. [4] 
The line profile of C ls O, which is assumed to be optically 
thin, looks roughly Gaussian. The only two exceptions are 
L 1455 Al and IRAS 03282, which show more complex line 
shapes, indicating that these sources may contain several veloc- 
ity components. The C I8 line width (FWHM) is quite narrow, 
typically <1 km/s. Such narrow line widths indicate that the 
bulk of the molecular gas is cold, and turbulence is low. The 
column densities we calculated for C I8 (see below), indicate 
that the C 18 1-0 line is optically thin in these objects. 

The line shape of HCO + 3-2 deviates clearly from a 
Gaussian profile. In 12 objects double-peaked profiles appear, 
with minima at the position of the optically thin C 18 emission, 
a clear indication for self-absorption. Asymmetries in the line 
shape indicate a systematic velocity pattern, which is probably 
caused by infall motions. 

Gregersen et al. (12000) observed several young protostars 
in HCO + 3-2 using the Caltech Submillimeter Observatory 
(CSO). They found line shapes that show similar asymmetries 
to those we detect (see section 14. 5t , but the intensities of the 
lines they observed are significantly lower. This discrepancy 
is most likely due to the larger beam size (32.5") of the CSO 
at 268 GHz. For three sources common to the Gregersen et 
al. (|2000) study and ours, we determined beam filling factors 
based on the HCO+ 3-2 lines. We used a CSO beam efficiency 
of 69.8% at 269 GHz to convert their T^ to Tmb- Assuming 
that HCO + 3-2 is emitted from a small region at the centre of 
the core, as suggested by our models (see section l5~4l i. we es- 
timated the size of this region by comparing the T MB ratio to 
the ratio of the beam size of the two observations. The derived 
filling factors for the CSO observations were between 0.26 and 
0.36, yielding a diameter of the warm nucleus of approximately 
4000 AU. Contributions of the error beams to the IRAM obser- 
vations have been neglected here, but due to the small size of 
the HCO + emitting region, their contribution is < 1%. 

4. Analyses 

In this section, we calculate the N2D + /N2H + ratios and corre- 
late them with physical parameters of the protostars that might 
influence deuterium fractionation. Some of these parameters, 
such as dust temperature and bolometric luminosity, can be 
used to quantify the evolutionary stage of the protostar. Each 
parameter will be discussed in individual subsections. 

4.1. N 2 H + & N 2 D + volume densities, column densities 
and deuterium fractionation 

To derive the deuterium enhancement, one calculates the ratio 
of the column densities of a hydrogen-bearing species and its 
deuterated isotopologue. Many common C-bearing molecules 
(e.g. HCO + and H2CO) suffer from depletion in cold and dense 
environments (see Tafalla et al. 2006), and thus the deuterated 
species trace slightly warmer regions of the core, where the 
parent molecules are not depleted. To properly determine the 



Table 4. Observational parameters of the C ls O and HCO + 
lines. Integrated HCO + intensities from double-peaked spectra 
are indicated with an *. 



Source 


W [KKm/sJ 


Av [km/s] 


W [JsJOTl/Sj 






(.<- <-)) 


) 


L 1448 IRS 2 


3.23±0.05 


0.77+0.01 


13.0+0.3* 


L 1448 IRS 3 


5.75±0.07 


1.22+0.01 


29.0+0.3 


L 1448 C 


3.36±0.06 


1.11+0.01 


20.6+0.3* 


L 1455 Al 


1.65+0.06 


1.03+0.01 


8.8+0.3 


L 1455 A2 


0.90±0.05 


0.80+0.01 


6.5+0.3 


Per4S 


2.84±0.06 


0.81+0.01 


8.7+0.3 


Per4N 


l.30±0.04 


f\ A A . f\ f\ 1 

0.44+0.01 


5.3+0.3 


Per 5 


1.79+0.05 


0.69+0.01 


15.0+0.3 


Per 6 


1.52±0.05 


0.74+0.01 


6.1+0.3* 


IRAS 03282 


1.95+0.07 


1.47+0.01 


12.8+0.3* 


Barnard 1-B 


8.09+0.08 


1.59+0.01 


14.0+0.3* 


HH211 


5.19+0.06 


0.81+0.01 


17.4+0.3* 


Barnard 5 IRS 1 


3.99+0.06 


1.00+0.01 


9.3+0.3* 


L 1527 


4.12+0.06 


0.65+0.01 


13.0+0.3* 


L483 


5.12+0.06 


0.81+0.01 


22.8+0.3* 


SMM5 


8.66+0.08 


1.57+0.01 


10.7+0.3 


L723 


2.45+0.08 


1.29+0.01 


12.2+0.3* 


L673 A 


2.31+0.07 


0.92+0.01 


12.8+0.3 


B335 


2.34+0.04 


0.44+0.01 


3.6+0.3* 


L 1157 


1.57+0.05 


0.61+0.01 


13.4+0.3* 



deuterium fractionation of an entire core, one has to observe 
molecules that are less affected by depletion, such as N2H + and 
NH3. Another benefit of observing N-bearing molecules is that 
their emission lines split into hyperfine components due to the 
non-zero nuclear spin of nitrogen, enabling optical depths to be 
measured. 

To determine molecular column densities, we used the 
N2H + 1-0 and N2D + 2-1 lines. We chose these lines for three 
reasons. First, N2D + 2-1 is detected in most of the sources. 
Second, the 2-1 transition of N2D + is observed with the high- 
est S/N ratio. Third, the ratio of the N 2 D + 2-l/N 2 H + 1-0 lines 
permits a direct comparison with previous studies (Crapsi et 
al. 120051 Fontani et al. 2006). The column densities were de- 
rived (Tab. |5]l assuming a constant excitation temperature, as 
determined in section[3](CTEX-method, for details see the ap- 
pendix of Caselli et al. 2002b). The excitation temperature of 
N2H + can be determined in all objects but Per 6 (Tab. [3]), where 
the N2H + line is optically thin. For Per 6, we assumed an exci- 
tation temperature of 8 K, which is the mean excitation temper- 
ature of the sample. For N2D + , the situation is different, since 
the N2D + 2-1 line is optically thin in most objects. In the few 
sources where r could be determined, the excitation temper- 
atures of N2H + and N2D + were approximately the same. In 
all cases where the excitation temperature of N2D + was un- 
known, we assumed that it was the same as for N2H + . The 
estimated errors on the column densities are of the order of 
15% and 20% for N 2 H + and N 2 D + , respectively. As expected, 
the N2D + /N2H + ratios vary significantly from object to object 
(Tab.|5]l. The extreme values are < 0.029 for Per 6 and 0.27 for 
HH 211. For most objects, the ratio is below 0.1. 
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Determining T ex and the column densities as described 
above, we implicitly assume a unity beam filling factor for the 
N2H + and N2D + observations. Because dark cloud cores are 
typically larger than the beam size of our observations (e.g. 
Caselli et al. 2002a) and the fact that the kinetic temperature 
of these cores is only ~ 1.5 times the excitation temperature 
of N2H + indicates that the filling factor is indeed close to one, 
and certainly > 0.67. Using the lower limit for the filling factor 
would lower the N2H + column density by 40% and the N2D + 
column density by 55%. Therefore, the derived N2D + /N2H + 
ratios could be up to 25% lower. The general trend, however, 
which we will discuss in the next section, remains. 

The N2D + column density was also calculated using the 
(3-2) line, where detected, following the same procedure de- 
scribed above. In 13 out of 15 objects, the values for the col- 
umn density agree within the range of error. In L 1448 IRS 3 
and L 1157, the column densities determined using N2D + 3-2 
are a factor of two higher than those determined using the 2-1 
line. However, because of the lower S/N ratio of the N2D + 3-2 
line, the column densities determined using the J =2-1 transi- 
tion are more reliable. 

Table 5. N 2 H + and N 2 D + column densities. 



Source 


N(N 2 H + ) 


N(N 2 D + ) 


N(N 2 D + )/N(N 2 H + ) 




[cm -2 ] 


[cm" 2 ] 




L 1448 IRS 2 


3.1-10° 


1.3-10 12 


0.042±0.014 


L 1448 IRS 3 


4.8-10 13 


3.9-10 12 


0.081+0.010 


L 1448 C 


2.4 10 13 


2.5-10 12 


0.104±0.017 


L 1455 Al 


9.2-10 12 


6.310" 


0.068±0.043 


L 1455 A2 


1.4-10 13 


9.6-10" 


0.069±0.029 


Per4S 


7.9-10 12 


1.4-10 12 


0.177±0.032 


Per 4 N 


7.1-10 12 


7.8-10" 


0.110±0.057 


Per 5 


1.2-10 13 


1.0-10 12 


0.083±0.033 


Per 6 


6.1-10 12 


<1.810" 


< 0.029 


IRAS 03282 


1.6-10 13 


3.5-10 12 


0.219±0.026 


Barnard 1-B 


4.5-10 13 


8.1-10 12 


0.180±0.009 


HH211 


1.7-10 13 


4.6-10 12 


0.271 ±0.024 


Barnard 5 IRS 1 


1.6-10 13 


5.4-10" 


0.034+0.025 


L 1527 


7.2-10 12 


< 2.5-10" 


<0.034 


L483 


3.3-10 13 


1.9-10 12 


0.057±0.012 


SMM5 


2.9-10 13 


2.3-10 12 


0.079+0.019 


L723 


8.8-10 12 


4.0-10" 


0.045+0.040 


L673 A 


1.3-10 13 


4.3-10" 


0.033+0.031 


B335 


1.410 13 


6.0-10" 


0.043+0.029 


L 1157 


1.1-10 13 


5.0-10" 


0.045+0.036 



Estimated errors on the column densities are about 15% and 20% for 
N 2 H + and N 2 D + , respectively. 

The N2H + and N2D + column densities derived here are 
of the same order as those found by Crapsi et al. (120051 1 for 
prestellar cores (10 13 cm -2 and 10 12 crrT 2 forN2H + and N2D + , 
respectively). This result is somewhat unexpected, because the 
optical depth of N2H + is about three times lower in proto- 
stars than in prestellar cores. However, the higher excitation 
temperatures and broader line widths of N2H + in protostellar 
envelopes (over-)compensate for the lower optical depth. The 



range of values for the N2D + /N2H + ratio is also similar. The 
average ratio in our sample of Class protostars is 0.097, and 
20% of the objects have a ratio larger than 0.15. For starless 
cores, the average ratio is 0.106, and 18% of the sources have 
ratios larger than 0.15. However, the N2D + /N2H + ratio in both 
protostellar and prestellar cores is significantly larger than the 
ratio in high mass star-forming regions (always smaller than 
0.02, Fontani et al. 2006), although values as large as 0.1 have 
been measured with interferometric observations (Fontani et 
al. l2008l) . 

For three objects in our sample, L 1448 C, HH 211, 
and L 1527, the N2D + /N2H + ratio was also determined by 
Roberts & Millar ( 2007 ). In general, they find somewhat higher 
N 2 D + /N 2 H + ratios than we do. In L 1448 C, their ratio is 
1 .63 +0.4 times larger. In HH 211, their ratio is 1.14+0.18 times 
larger. For L 1527, we report an upper limit of 0.034, whereas 
Roberts & Millar found a ratio of 0.06. The fact that Roberts & 
Millar (120071) found similar N 2 D + /N 2 H + ratios at much lower 
spatial resolution indicates that the deuterium fractionation oc- 
curs, as expected, in the cold, extended envelope of the pro- 
tostellar core. If deuterium fractionation occurred only in the 
central region, one would expect beam dilution to yield lower 
N2D + /N2H + ratios at lower resolutions. 

We also estimated the volume densities of the cores. 
Because the N2D + lines show only a single velocity compo- 
nent, we used the ratios of the integrated intensities of these 
lines instead of the Tmb ratios. Thus we obtained a higher sig- 
nal to noise ratio. Further we assumed an isothermal cloud 
with a kinetic temperature determined from NH3 observa- 
tions (Tab. [3). The volume density was estimated by fitting 
the N2D + J=3-2/J=2-l ratio, using a radiative transfer code 
(Pagani 2007, see section |5~2| |. We derived values between 
5 ■ 10 6 cm' 3 and 10 7 cm 4 . Only IRAS 03282 shows a lower 
density (6 • 10 5 cm -3 ). However, temperature and density de- 
cline with distance from the protostar, so lower spatial reso- 
lution observations trace lower average temperatures, densi- 
ties and (assuming spherical symmetry) lower column densi- 
ties. Therefore, the N2D + 3-2 lines trace, on average, warmer 
and denser gas than the N2D + 2-1 lines, and the densities ob- 
tained from the ratios are upper limits. For six sources, we also 
determined densities using the N2D + J=2-l/J=l-0 ratios. We 
were not able to find a density which would reproduce all three 
lines simultaneously. The density values derived from the J=2- 
l/J=l-0 line ratio are about a factor 10 lower than those de- 
termined from the J=3-2/J=2-l ratio. This difference is again 
most likely arises because of the lower spatial resolution of the 
lower frequency observations. 

4.2. Dust Temperature 

In general, Toust is expected to rise as a protostar evolves, 
but there is no exact model of dust temperature as a func- 
tion of time. As a result, connecting dust temperature to spe- 
cific protostellar ages depends on which model is used (see 
Froebrich 2005 and references therein). 

All of our sources, with the exceptions of Per 4 S and 
SMM 5, have infrared counterparts identified in the IRAS point 
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source or faint source catalogues. In addition, millimetre and 
submillimetre continuum data are available for many of the 
sources (see Tab. [6j. We fitted these data to derive values for 
dust temperature (Tpust) and emissivity (J3) using equation (1) 
from Kramer et al. (20Q3J: 

1 



F v oc v 



3+j8 



e hv/kT Dusl _ J 

where F v is the continuum flux at frequency v. For 15 out of 
16 sources for which more than two data points are available, a 
value of /3- 1 yielded the best fit. Only for Barnard 5 IRS 1 does 
an index of 0=0 (i.e. black-body radiation) yield a better re- 
sult. In the three cases where data at only two frequencies were 
available, we assumed /3- 1 . This value is lower than the typical 
emissivity index of f3 — 2 for molecular clouds (e.g. Ward- 
Thompson et al. [2002, Visser et al. 1998). We propose that 
the lower emissivity indices found in the Class sources are 
a result of strong temperature gradients and not inherent differ- 
ences in the dust grains themselves. In fact, the cold envelopes 
are still bright at millimetre wavelengths, whereas practically 
all radiation emitted at A < 100/im comes from the hot core. 
Therefore, the mass of dust radiating at longer wavelengths is 
much greater than the mass of dust radiating at shorter wave- 
lengths, and the spectrum from the entire core gets flattened. 
We modelled the continuum emission using an emissivity index 
of two (see section[5]l, but got as index for the whole modelled 
cloud a /3-value of one. 




Fig. 3. Deuterium fractionation (N2D + /N2H + column density 
ratio) versus dust temperature. The largest N2D + /N2H + ratios 
are seen at the lowest temperatures. As Toust increases, the deu- 
terium fractionation declines. The objects marked with trian- 
gles are located in the Perseus cloud. 

The dust temperatures derived from fitting the continuum 
data range from approximately 20 K up to 50 K (see Tab. [6j. 
These temperatures are clearly higher than the excitation tem- 
peratures of N2H + , as well as the kinetic gas temperatures 
(Tab. [3}. The difference is likely due to the fact that the dust 
emission from warmer regions is stronger than from cold dust, 
and thus the obtained average temperature is higher than the 
gas temperature traced by N2H + . 

In Fig. [3] the N2D + /N2H + column density ratio is plotted 
versus the dust temperature. All objects with N2D + /N2H + ra- 



tio >0.15 have Tn U st < 25 K. At higher dust temperatures, the 
deuterium fractionation declines, reaching an average value of 
approximately 0.05 at 30 K. 

4.3. CO freeze-out 

Because CO is one of the main destroyers of H3 and deuter- 
ated isotopologues, the progenitors of N2D + , the N2D + /N2H + 
ratio should be sensitive to CO depletion (e.g. Dalgarno & 
Lepp [19841 Roberts & Millar EOOOl The CO depletion factor, 
f D , is defined as 

— Xcan/^obs 

where x can is the canonical ISM abundance ratio of CO. In this 
work, C ls O has been observed, so that x can = [C 18 0]/[H 2 ] = 
1.7 • 10~ 7 (Frerking et al. I1982K The total number of hydro- 
gen molecules within the beam is determined by the ~ 1 mm 
continuum emission (1.3 mm, 1.1 mm, or 850 /mi), using the 
equation 

N F v ■ d 2 1 

K v ■ B v (TDust) H 

(Terebey et al. 119931 1. where F v is the continuum flux, d is the 
distance to the object, B v (To U st) is the Planck function, To U st is 
the dust temperature determined in the previous subsection, and 
H is the mass of an H2 molecule. For the opacity per unit mass, 
k v , we used a value of 0.005 cm 2 /g at a wavelength of 1.2 mm 
(Andre et al. 1996). To calculate k v at other wavelengths, we 
express the relation between opacity and wavelength as a power 
law 



We used ft = 1 (zero for Barnard 5 IRS1), as suggested by 
the fits to the continuum data (see section l4~2"b . The values for 
k make the largest contributions to the uncertainties in the H2 
column densities. Some authors use values for k twice as high 
as ours (e.g. Froebrich|2005ll. 

To compare the ~ 1 mm continuum data with our C ls O 1- 
observations, we have to estimate the continuum flux den- 
sity within a 22" beam. The fact that the angular diameters of 
the protostellar envelopes are always larger than the beam is 
an important issue in this conversion. To compare data taken 
with different beam sizes, we follow the method of Terebey et 
al. ([1993 ). They showed that the flux of a dark cloud with a den- 
sity distribution p = r~ p and a temperature distribution T = r~ q 
can be described by: 

F v oc d Mp+q) • 6» 3 - (p+q) 

where d is the distance and 6 is the beam size. In their sample 
of mostly Class protostars, they determined 1 .7 to be an ap- 
propriate value for (p+q). We used this value to scale the flux 
from the continuum observations to the C 18 beam size. 

We determined the C ls O column density by using the 
CTEX-method and by assuming that C 18 emission is opti- 
cally thin. Since the critical density of C 18 is 60 times lower 
than that for N 2 H + , we expect that T ex (C ls O) will be higher 
than T ex (N2H + ) and comparable to the kinetic gas temperature 
or even slightly higher, due to freeze-out of CO in the cold 
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Table 6. Continuum observations and derived properties (dust temperature, C 18 column density, CO depletion factor) for 
sources in our survey. The errors in the dust temperatures are ~ 2 K, and the C 18 column densities are accurate to 10 %. 



Source 


IRAS 


F(60yum) 


F( 100pm) 


Wavelength of additional data 


TDust 


N(C ls O) 


f D (CO) 






[Jy] 


[Jy] 


w 


[K] 


r -?n 

[cm 2 ] 




L 1448 IRS 2 


03222+3034 


14.3 


<169 


450, 850, 1100 u 


27 


6.07-10 15 


2.6 ± 0.6 


L 1448 IRS 3 


03225+3034 


25.1 


169 


1100 3 


27 


1.08-10 16 


4.4 ± 1.0 


L 1448 C 


F03226+3033 


53.4 


355 


450, 850, 1100 2 - 3 


32 


7.18.19-10 15 


2.8 ± 0.6 


L 1455 Al 


03245+3002 


47.1 


93.6 


1100, 450, 850 3 ' 5 


37 


3.53-10 15 


1.5 ±0.4 


L 1455 A2 


F03247+3001 


66.4 


206 


- 


30 


1.83-10 15 


- 


Per4S 


- 


- 


- 


1100 3 


- 


3.18-10 15 


- 


Per4N 


03262+3123 


1.07 


<12.4 


1100 3 


23 


2.18-10 15 


2.4 + 0.11 


Per 5 


03267+3128 


1.88 


<10.4 


1100 3 


26 


3.28-10 15 


1.6 ±0.4 


Per 6 


03271+3013 


7.53 


8.19 




45 


4.28-10 15 




IRAS 03282 


03282+3035 


2.33 


13.9 


450, 850, 1100 2 - 3 


23 


3.28-10 15 


3.9 ±0.9 


Barnard 1-B 


03301+3057 


<6.18 


35.5 


1100, 850 3 ' 4 


18 


1.1610 16 


2.9 ±0.7 


HH211 


03407+3152 


<2.06 


<30.8 


1100, 450, 850 3 - 5 


<21 


8.21-10 15 


2.9 ±0.7 


Barnard 5 IRS 1 


03445+3242 


15.5 


15.4 


1100 3 


50 


1.23-10 16 


0.3 ± 0.09 


L 1527 


04368+2557 


17.8 


73.3 


450, 850, 1300 2 


27 


7.75-10 15 


2.3 ±0.5 


L483 


18148-0440 


89.1 


166 


450, 850 2 


34 


1.15-10 16 


2.1 ±0.5 


SMM5 












9.70- 10 15 




L723 


19156+1906 


6.93 


20.7 


450, 850 2 


25 


4.36-10 15 


2.9 ±0.7 


L673 A 


19180+1114 


2.55 


10.6 




28 


4.47- 10 15 




B 335 


19345+0727 


8.3 


42 


450, 850, 1300 2 


27 


4.41-10 15 


3.3 ±0.8 


L 1157 


20386+6751 


10.9 


43.5 


1300, 450, 850 2 ' 5 


30 


3.19-10 15 


4.2 ± 1.0 



'O'Linger et al. ( 1999 1, 2 Shirley et al. ( 2000), 3 Enoch et al. (2006), 4 Hatchell et al. {2005), 5 Froebrich et al. ( [20031 . 



parts of the cloud. Since the beam sizes of the NH3 observa- 
tions which were used to derive T^n vary between 40" and 88", 
we adopted a constant excitation temperature of 15 K for all 
sources. This assumption introduces some error (< 15%) in the 
estimation of the depletion factor. The total number of C ls O 
molecules is then calculated from 

Wc>*o = N c is • d 2 • Q 

where Mc'«o is the total number of C ls O molecules within the 
beam, N c is is the C 18 column density, d is the distance to the 
protostar, and Q is the beam solid angle. 

The derived depletion factors lie between 1 and 4 (Tab. [6]). 
A value less than 1 is found toward Barnard 5 IRS 1 and this 
may be an indication that the C 18 canonical abundance is 
larger in this source. Indeed, the 13 CO abundances measured 
toward Barnard 5 by Pineda et al. (120081) are a factor of two 
larger than those measured by Frerking et al. (1982) toward 
Taurus (CO abundance variations of about a factor of 2 are well 
known in star forming regions, e.g. Lacy et al. |1994| ). However, 
given that CO abundance measurements are not individually 
available for all the sources, the C 18 canonical abundance has 
been kept the same for the whole sample. 

In the central panel of Fig. |H we show the correlation be- 
tween CO freeze-out and deuterium fractionation in our sam- 
ple of protostars. In the top panel, we show the same corre- 
lation found in low-mass prestellar cores (Crapsi et al. [2005), 
and in the bottom panel, the correlation for massive prestellar 
cores (Fontani et al. 120061 ) . The correlation determined from 
our sample looks very similar to the correlation found for low- 
mass prestellar cores. For small depletion factors (fp(CO)< 3), 



the deuterium fractionation is low as well (< 0.1). For depletion 
factors greater than three, the N2D + /N2H + ratio rises quickly 
until it reaches 0.25. The only difference between the protostel- 
lar and prestellar cores is that in the prestellar cores, the deu- 
terium fractionation stays low until the depletion factor is ~ 10. 
This behavior reflects the fact that CO is highly depleted in the 
center of prestellar cores, whereas it is not depleted in the cen- 
ters of Class sources. In high mass prestellar cores, there is 
no dependence of the deuterium fractionation on fo(CO). Since 
the high mass prestellar cores in the Fontani et al. 2006 sample 
are at larger distances (a few kpc), these observations might be 
affected by non-depleted, non-deuterated gas along the line of 
sight. 

The sources that are not located in the Perseus cloud show 
lower deuterium fractionation than Perseus sources with the 
same CO depletion factor. One explanation for this difference 
is that the Perseus sources are embedded in a dense environ- 
ment (ne 2 ~ 10 4 crrT 3 , Kwon et al. 2006), whereas many of 
the other cores are isolated. Since CO is not expected to be 
frozen out in the embedding material, the CO depletion fac- 
tors of the Perseus cores could be underestimated. However, 
systematic errors arising from the different spatial resolutions 
and different wavelengths of the continuum observations can- 
not be ruled out. The Perseus cores were observed at 1 100 ym 
with a resolution of 3 1 " (Enoch 2006); the other cores were ob- 
served at 1 300 \im with a resolution of 40" (Shirley et al. 120001 ). 
The Perseus source which does not appear to agree with the 
correlation is L 1448 IRS 3. It has a CO depletion factor of 
4.4 but a low deuterium fractionation (0.081). This source is 
known to consist of multiple cores with several strong outflows 
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Fig. 4. N2D + /N2H + column density ratio versus CO depletion 
factor for prestellar cores (upper panel, Crapsi et al. 2005), 
Class protostars (middle panel, this work) and massive 
prestellar cores (lower panel, Fontani et al. |2006l l. The proto- 
stellar objects marked with triangles are located in the Perseus 
cloud. Among the Perseus objects, fo(CO) and the deuterium 
fraction appear correlated. The solid curve in the middle panel 
is the result of fitting the sub-sample of Perseus cores with a 
formula adopted from Caselli et al. (1998 , see text). 



(e.g. Volgenau et al. 120061 Looney et al. 2000), but it is not 
clear how this activity can decrease the deuterium fractionation 
while maintaining a large amount of CO freeze-out. 

The correlation among the sub-sample of Perseus sources 
is quite good, whereas the non-Perseus sources do not show 
this functional dependence. We made a least-squares fit of the 
Perseus data to the equation 

N 2 D + /N 2 H + = , 

(c 2 +k co -x can (CO)/f D (CO))' 

which we adopted from expression (1) in Caselli et al. {1998 ). 
In this equation, Ci = 1/3 • kf • x(HD), (where kf is the for- 
mation rate of H 2 D + and x(HD) is the HD abundance), c 2 is 
the electron abundance (x(e)) times the dissociative recombi- 
nation rate of H 2 D + , k co is the destruction rate of H 2 D + by CO, 
and x can (CO) is the canonical CO abundance. The equation ne- 
glects the effect of multi-deuterated species and also assumes 



that dissociative recombination is a significantly faster process 
than recombination onto dust grains. A more comprehensive 
chemical model will be used in Sect. [5] The result of the least- 
squares fit is indicated by the solid curve in Fig. [4] The cor- 
relation coefficient is 0.71, indicating that the data agree fairly 
well with the expected theoretical fit. Using the reaction rates 
given by Caselli et al. (1998) and assuming a temperature of 
10 K, we determined the HD abundance and the electron abun- 
dance. The derived an HD abundance of 3 • 10~ 5 is, as expected, 
quite close to the cosmic D/H ratio. The electron abundance we 
find is ~ 7 • 10~ 9 , which is at the low end of the values found 
in dense cloud cores. However, due to the uncertainties in the 
temperatures, densities, and the determination of fo(CO) (see 
below), errors in the HD and electron abundances determined 
in this way may be as large as an order of magnitude. 

The errors in the depletion factor are of the order of 20 %. 
These errors do not include the uncertainties in the excitation 
temperatures, dust temperatures and /c's, which are the main 
sources of systematic errors in this analysis. A different k value, 
however, would just rescale fo, but it would not change the 
correlation itself. Such a rescaling would also be the effect of 
a different dust temperature, because fp(CO) was determined 
from continuum measurements around 1 mm, which is a wave- 
length clearly in the Rayleigh-Jeans regime. Rejecting the as- 
sumption of a constant average excitation temperature in all 
objects would change the correlation significantly. 

If T ex (C ls O) is lower than 15 K, we Mnoferestimate 
the depletion factor, if T ex (C ls O) is greater than 15 K we 
overestimate it. But, if we assume that T ex (C ls O) is related to 
the dust temperature, i.e. objects with higher Toust have also 
higher C ls O excitation temperatures, and taking the correla- 
tion between the N 2 D + /N 2 H + ratio found in the previous sec- 
tion into account, the high deuterium fractionation cores would 
may have greater fo(CO) values, and the low deuterium frac- 
tionation cores would have lower fo(CO) values. Therefore the 
correlation would spread out further along the depletion factor 
axis. 



4.4. Bolometric Luminosity 

The SED of a young stellar object changes as the central pro- 
tostar heats up the surrounding gas and dust. Several mea- 
sures have been developed to determine how far a proto- 
star has evolved, including the bolometric temperature (Myers 
& Ladd 1993 ), the ratio of the submillimetre to bolomet- 
ric luminosity (Andre et al. 119931) . and the correlation be- 
tween luminosity (L) and visual extinction (A v ) (Adams Tl 9901 ). 
Another measure of the evolution of a protostar, developed by 
Saraceno et al. ( 1 19961 ), is the ratio of Lbol to FA, the distance- 
normalised flux at 1.3 mm. As the protostar warms up, Lbol 
should rise quickly, whereas F^, which is a measure for the 
envelope mass and less dependent on the temperature, should 
stay constant or decrease slightly. Hence, the ratio Lbol/F°3 
increases as the protostar evolves. 

To determine F^ values, we use fluxes at A = 
1.3 mm (5 sources), 1.1 mm (4 sources) and 850 /mi (2 sources). 
The continuum fluxes at 1 . 1 mm and 850 //m were corrected 
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using, as before, (3=1 (or, in the case of Barnard 5 IRS 1, 
ft = 0). For the distance normalisation and the conversion of 
the different beam sizes, we used the equation from Terebey et 
al. ([1993) described in the previous section, i.e. 

F 220 P c _ F obs. . (d/2 20pc) - 7 • (6>/31"r L \ 

where F° b 3 s is the continuum flux, d is the distance to the ob- 
ject, and 6 is the beam size. We normalised all observations to a 
distance of 220 pc and a beam size (FWHM) of 31". From the 
errors given for the continuum fluxes and bolometric tempera- 
tures, we estimate an error of 30 % in the calculated Lbol/F^ 3 
ratios. The resulting normalised fluxes and the bolometric lu- 
minosities for the cores are listed in Tab. [7] 

Table 7. Bolometric luminosity and distance-normalised flux 
of Class protostars in our sample. 



Source 


Lbol [L©] 


R OT pc 
M.3 


T fC I2U P c 


L 1448 IRS 2 


6.0' 


0.80 


7.5 


L 1448 C 


8.3' 


0.53 


15.6 


L 1455 Al 


9.7 1 


0.19 


50.2 


IRAS 03282 


1.2 2 


0.32 


3.8 


Barnard 1-B 


2.8 


0.74 


3.8 


HH211 


3.6 1 


0.61 


5.9 


L 1527 


1.9 1 


0.38 


5.0 


L483 


13.0 2 


0.70 


18.5 


SMM5 


3.6 1 






L723 


3.3 2 


0.45 


7.4 


B335 


3.1 2 


0.45 


6.9 


L 1157 


5.8 2 


0.68 


8.6 



'Froebrich (2005), 2 Shirley et al. f2000t . 



all have low Lbql/Fj 3 ratios and are presumably less-evolved 
protostars. As the Lbol/F° 3 ratio increases, the deuterium frac- 
tionation drops quickly. If one considers just the Perseus cores, 
then the fall-off of the deuterium fractionation is less steep. The 
L B ol/F°3 ratios may be lower in the non-Perseus sources be- 
cause systematic errors in normalisation (see section l4~3l lead 
to overestimates of the continuum fluxes. L 1455 Al has by 
far the highest Lbol/F° 3 ratio, and by this measure should 
be the most evolved protostar in the figure. The assertion that 
L1455 Al is more evolved that other protostars is reinforced by 
the relatively low CO depletion factor and high dust tempera- 
ture. L 1527 is a notable outlier in the sample. Its Lbol/F° 3 
ratio is one of the lowest in our sample (5.0), but it has an 
N2D + /N2H + ratio < 0.034 (more discussion on this source is 
presented in 14. 1\ . 

4.5. Kinematic of the Gas 

So far, we have correlated the N2D + /N2H + ratios with physical 
parameters, which are all related to temperature, and thus to the 
evolutionary stage of the protostar. In this section, we inves- 
tigate the relation of deuterium fractionation with parameters 
linked to the kinematics of the gas. Infall, for example, heats up 
the gas and is expected to influence the gas-phase chemistry. 

First, we compare the N2DVN2FP ratios with the N2D + 2-1 
line widths (Tab.0. Processes such as thermal broadening, tur- 
bulent motions, and systematic motions, e.g. infall, contribute 
all to the full width at half maximum (FWHM) of optically 
thin emission lines. Thermal broadening plays a relatively mi- 
nor role, since the thermal line width (FWHM) of N 2 D + at 10 K 
is only 0.12 km/s (0.17 km/s at 20 K). 
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Fig. 5. N 2 D + /N 2 H + versus L BO l/F°° 3 . All objects with high 
deuterium fractionation have low Lbol/F° 3 ratios and are 
likely to be young objects. L 1527, which has one of the 



lowest Lbol/F° 3 ratios, also has a quite upper limit on the 
N 2 D + /N 2 H + ratio. 



The correlation between 



Lbol/F°- 



3 and N 2 D + /N 2 H + is 
shown in Fig. [5] The objects with the highest deuterium frac- 
tionation (> 0.15) (IRAS 03282, Barnard 1-B and HH 211) 
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Fig. 6. FWHM of N 2 D + 2-1 versus N2D + /N 2 H + . Triangles are 
Perseus cores. 

In Fig. [6] the N2D + /N2H + column density ratio is plotted 
versus the FWHM of N 2 D + (2-1). Most of the protostars with 
high N2D + /N2H + ratios have narrow N2D + lines. The two ex- 
ceptions are Barnard 1-B (Av = 0.78 km/s) and L 1448 C 
(Av = 0.72 km/s), which show the broadest N2D + lines in the 
sample. However, all Perseus sources with Av > 0.5 km/s, in- 
cluding Barnard 1-B and L 1448 C, have multiple embedded 
sources (Volgenau 2004, Matthews & Wilson l2502b . The rela- 
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tive velocity differences between the individual protostars leads 
to additional line broadening. 

To be able to trace only the systematic motion of the gas, 
e.g. infall motion, we used the asymmetry parameter Sv, which 
was defined by Mardones et al. d 19971 ) as 



Sv 



v thick ~ v thin 
AVthin 



Vthick is the velocity of an optically thick line, v t h; n is the ve- 
locity of an optically thin line, and Av t hi n is its line width. 
As shown by Mardones et al., infall yields a negative value 
for <Sv. If infall dominates the kinematics of the protostar, then 
an absorbing foregroud component will appear redshifted with 
respect to the core, causing an optically thick line to appear 
blueshifted. We used HCO + 3-2 as the optically thick line and 
N2H + 1-0 as the optically thin line. The velocity of the optically 
thin N2H + line is derived using the HFS method of CLASS (see 
section[3]l. We also compare the resulting velocities with those 
derived by fitting a single Gaussian to the isolated lot — > 1 12 
component, following Mardones et al 1997 The velocities cal- 
culated with the two different methods are not always coinci- 
dent, but the differences are within the velocity resolution of 
0.06 km/s. Because a fit to seven lines has a lower error than 
a fit to only one hyperfine component we use in the following 
the Vlsr obtained by the HFS method. For v t h; c k, we use the 
peak velocity of the HCO + 3-2 line. In cases where the HCO + 
lines are double-peaked, we follow Mardones et al. ( 1997 ) and 
use the velocities of the stronger of the two peaks (although 
in some spectra the stronger peak was ambiguous). An alter- 
native asymmetry parameter is the ratio of the red peak inten- 
sity to the blue peak intensity (R/B) of optically thick, double 
peaked lines. In such cases, R/B < 1 indicates infall. Since 
more than 50 % of the observed HCO + lines show double- 
peaked line shapes, we also calculate these ratios. The results 
for both asymmetry parameters are listed in Tab. [8] 

In our sample, we find eight sources with a clear indication 
of infall (i.e. <5v < -0.25), and seven sources with redshifted, 
optically thick lines (i.e. Sv > 0.25). For five sources, there is 
no clear indication of systematic motion (-0.25 < Sv < 0.25). 
Like previous studies by Mardones et al. {T997 ) and Gregersen 
et al. ( 120001) . we find more blueshifted than redshifted objects. 
The blue excess is defined by Mardones et al. (|T997 ) as 

Nblue - N re d 



blue excess = 



N total 



where N D i ue is the number of objects with Sv < -0.25 and N re d 
is the number of objects with Sv > 0.25. In our sample, the 
blue excess is 0.05 + 0.2. If we exclude the objects with an 
ambiguous Sv due to nearly equal peak intensities of the double 
peaked HCO + line, we derive a blue excess of 0.12 + 0.2. Both 
values are slightly lower, but still agree within the errors with 
0.28 and 0.25, the values found by Gregersen et al. (2000) and 
Mardones et al. ( 1997), respectively. 

Sv and, less clearly, R/B are correlated with the deuterium 
fractionation, and only Per 5 and, once again L 1527, which 
is a peculiar source (see section 14.71 ). do not follow the gen- 
eral trend and show indications of strong infall motion. The 



Table 8. Velocities of C ls O and HCO + lines and asymmetry 
parameters. Values for Sv marked with * are ambiguous due to 
nearly equal peak intensities of the double-peaked HCO + line, 
i.e. the R/B ratio is close to one. 



Source 


v [km/s] 
N 2 H+ 


v [km/s] 
HCO + 


<5v 


R/B 


L 1448 IRS 2 


4.12 


4.53 


0.82* 


1.28 


L 1448 IRS 3 


4.62 


5.40 


0.84 


- 


L 1448 C 


4.95 


5.45 


0.66 


2.21 


L 1455 Al 


5.10 


5.12 


0.03 


- 


L 1455 A2 


4.88 


4.94 


0.17 


- 


Per4S 


7.52 


7.25 


-0.92 




Per4N 


7.57 


7.44 


-0.40 




Per 5 


8.15 


7.65 


-1.31 




Per 6 


5.94 


6.02 


0.05 


1.82 


IRAS 03282 


6.06 


6.62 


-0.91 


0.57 


Barnard 1-B 


6.62 


6.01 


-0.74 


0.45 


HH211 


9.12 


9.27 


0.37* 


1.47 


Barnard 5 IRS 1 


10.27 


10.29 


0.05 


3.88 


L 1527 


5.90 


5.38 


-1.60 




L483 


5.40 


5.87 


1.06 


3.18 


SMM5 


8.12 


8.18 


-0.33 




L723 


11.11 


11.43 


0.46 


1.61 


L673 A 


6.87 


6.97 


0.19 




B 335 


8.34 


7.90 


-1.10* 


0.833 


L 1157 


2.67 


3.16 


0.75 


1.54 




Fig. 7. Asymmetry parameter Sv versus N2D + /N2H + ratio. 
Only the sources with an unambiguously determined asymme- 
try parameter are shown. All sources with enhanced deuterium 
fractionation are located at moderate Sv of about -0.4, and is 
increasing with decreasing deuterium fractionation. 



higher the N2DVN2I-F ratio is, the more blueshifted the opti- 
cally thick emission lines are. Such a result was somehow ex- 
pected, since outflow activity rises with the evolution of the 
protostar (Richer et al. 2000). The results are also in good 
agreement with the results of Mardones et al. (1997). Using 
observations of CS, C 34 S and H2CO lines, they found that 
the blue excess of Class protostars is clearly enhanced com- 
pared to Class I objects. However, using observations of HCO + , 
H 13 CO + and N 2 H + lines, Gregersen et al. (2000) did not find 
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such trend. They speculated that infall asymmetry disappears 
later in HCO + . Because our sample contains no true Class I 
objects, we cannot directly compare their results with ours. 

We also calculated 5v using the C ls O 1-0 transition as the 
optically thin line. Contrary to N2H + , which traces the cold, 
CO depleted envelope, C ls O traces the already warmed up cen- 
tral part of the protostar. Thus it traces the region, where the 
HCO + 3-2 line stems from, too. But due to the lower criti- 
cal density of the C ls O 1-0 transition, the low density out- 
skirts of the protostellar envelope might contribute significantly 
to the C ls O emission, which possibly dilute any correlation. 
Qualitatively, the correlation between 5v and the deuterium 
fractionation looks similar, no matter which optically thin line 
is used. However, for objects with an N2D + /N2H + ratio above 
0.75, 5v determined with C 18 yields to less blue-shifted val- 
ues. A reason for this difference might be that in these objects, 
which also show the highest CO depletion factors, the contribu- 
tions from the environmental material to the C 18 emission is 
expected to be biggest. The bigger contribution of the environ- 
mental material is also indicated by the fact that for objects with 
Svn 2 h+ lower than -0.5, Av(C ls O) is on average 2.4 times larger 
than Av(N2H + ), whereas for the remaining sources Av(C ls O) 
is only 1 .6 times larger. More turbulent material appears to sur- 
round infalling protostellar envelopes. 

4.6. Comparison of the deuterium fractionation of 
different molecules 

The deuterium fractionation of N2H + reflects the present time 
temperature and density better than most other molecules, 
since it does not significantly deplete (see e.g. Caselli et 
al. I2002bl Bergin et al. [20021 Schnee et al. l2"0"07j l and the 
chemical reaction for N2D + formation and destruction are quite 
fast and simple. Memories of cold phases in the past of the 
dark cloud core are therefore not expected to be seen in the 
observed N2D + /N2H + ratio. In this subsection, we compare 
the N 2 D + /N 2 H + ratio with the NH 2 D/NH 3 (Hatchell 125531 . 
DCO + /HCO + ratios (J0rgensen et al. 2004) and the deuterium 
fractionation of H 2 CO and CH 3 OH (Parise et al. 12006b for the 
sources in our sample where these data are available (Tab. [9]). 

Table 9. Deuterium fractionation of several molecular species. 



Source 


N 2 D T l 


NH 2 D 2 


DCO+ 3 


D 2 CO 4 


N,H+ 


NH, 


HCO+ 


HbCO 


L1448 IRS 3 


0.081 


0.17 






L1448 C 


0.104 


0.20 


0.011 


0.24 


IRAS 03282 


0.219 


0.22 






HH211 


0.271 


0.33 






Barnard 5 


0.034 


0.18 






L 1527 


< 0.034 




0.048 


0.44 


L483 


0.057 




0.006 




L723 


0.045 




0.004 




L 1157 


0.045 




0.016 


<0.08 



'This work, 2 Hatchell (2003), 3 j0rgensen et al. l [2"004l , 4 Parise et 
al. 120061 . 



The deuterium fractionation of HCO + is, in all objects 
but one (L 1527), significantly lower than the fractionation 
of N2H + . This difference can be explained by the fact that 
CO, the parent species of HCO + , is also the main destroyer 
of H 2 D + and thus the DCO + /HCO + ratio traces only regions 
where the deuterium fractionation in N2H + is expected to be 
lower. Furthermore, J0rgensen et al. (2004) determined the 
DCO + /HCO + ratio using the DCO + J=3-2 line, which has a 
critical density of 1.4 • 10 6 cm 3 . Therefore this DCO + obser- 
vations traces mostly the inner, high density part, which might 
be already warmed up by the embedded protostellar core, thus 
lowering down the DCO + /HCO + ratio. The situation with NH3 
is different. Like N2H + , NH3 is formed from N2 and does not 
seem to freeze out (Tafalla et al l2002l . Thus, N 2 D + and NH 2 D 
emission should originate from the same region. As a conse- 
quence, the deuterium fractionation of N2H + and NH3 should 
be comparable. 




Fig. 8. Deuterium fractionation of NH3 (red triangles) and 
N2I-P (black circles) versus dust temperature. 

As shown in Fig.[8j the coldest objects, i.e. the youngest ob- 
jects in our sample (HFf 21 1 & IRAS 03282), have a deuterium 
fractionation of NH3 indeed comparable to the fractionation of 
N2H + . But the warmer the dust gets, the more the ratios di- 
verge. The NH2D/NH3 ratio is about 0.2 at T Dust = 50 K, where 
the N2DV N2PF ratio is five times lower. The differing trends 
in the two ratios can be explained as follows. Both N2D + and 
NH2D are formed directly or indirectly via H2D + or its multi- 
ply deuterated forms, i.e. the formation mechanism is only ef- 
fective in cold, CO-depleted environments. In regions were CO 
desorbs from the dust due to the warming up of the protostar, 
N2D + and N2H + disappear quickly from the gas phase, because 
these species are directly destroyed by CO, forming HCO + 
and DCO + , respectively. Thus, with increasing temperature, the 
N2D + /N2H + ratio traces material at larger and lager radii (i.e. 
at lower and lower densities). But deuterium fractionation is 
related to the density, because of the lower degree of CO de- 
pletion at lower densities (see Sect. 15.41 and Fig.fTTI). Contrary 
to N2H + and N2D + , which are pure gas-phase species, ammo- 
nia is also found on dust grain mantles, most likely formed by 
grain-surface reactions (Boogert et al. 120081 ). The solid NH3, 
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which desorbs form the dust grains as soon as the environmen- 
tal material warms up, likely has a large deuterium fractiona- 
tion due to the fossil record of grain surface chemistry and gas 
phase depletion at low temperatures (e.g. Aikawa et al.|2008). 
Thus, the observed NH2D/NH3 ratio declines more slowly than 
the deuterium fractionation of pure gas phase species. 

Deuterated formaldehyde (HDCO and D2CO) are very 
abundant, too. The HDCO/H2CO ratio as well as the 
D 2 CO/H 2 CO ratio are of the order of 0.1 (Parise et al. [2006t . 
In their work, Parise et al. also determined the deuterium frac- 
tionation of methanol, which exceeds even the deuterium frac- 
tionation of formaldehyde (e.g. CH 2 DOH/CH 3 OH ~ 0.5). 
For both molecules, no obvious correlation between the deu- 
terium fractionation and the core evolutionary stage could be 
found. Because mainly grain surface chemistry is involved in 
the deuteration of methanol and formaldehyde, the observed 
deuterium fractionation is likely a memory of the previous, 
prestellar phase (Ceccarelli et al. 2001), similar (although more 
extreme) case as found for ammonia. 

4.7. Evolutionary sequence of the sample 

The youngest objects in our sample are Barnard 1-B, 
IRAS 03282 and HH 211. All three evolutionary tracers indi- 
cate those sources to be among the youngest. The depletion fac- 
tor, f D (CO), indicates that L1448 IRS3 and L 1157 are young 
objects, too, but their Tqust suggest that these sources are, in- 
stead, more evolved. A possible solution to this discrepancy lies 
in acknowledging several assumptions we made in determining 
the CO depletion factor (e.g. in converting the beam sizes of 
the continuum and the C ls O observations). Although, overall, 
these assumptions seem to be justified, they might not be valid 
for all sources and thus the fo (CO) values may be less accurate. 
The three protostars we identified as the youngest also have the 
highest N 2 D + /N 2 H + ratios (0.18, 0.22 and 0.27 in Barnard 1- 
B, IRAS 03282 and HH 211, respectively). By this measure, 
Per 4 S, with a deuterium fractionation of 0. 177, should also be 
very young, but no complementary data are available to con- 
firm this. 

The objects at a moderate evolutionary stage are those with 
a deuterium fractionation around 0.1 (L 1448 IRS 3, PER 4 N, 
PER 5). The exception is L 1448 C (N 2 D + /N 2 H + = 0.10), 
which seems to be more evolved than the others. However, 
among these intermediately-evolved sources, there are also 
some with deuterium fractionations of ~ 0.05 (L 723, B 335, 
L 1448 IRS 2). For SMM 5 (N 2 D + /N 2 H + = 0.079), we can't 
determine the evolutionary stage independently, but, because 
of its N2D + /N2H + ratio, we think that it belongs to the group of 
intermediately-evolved objects. The N2D + /N2H + ratios of the 
most evolved protostars are, with the exception of L 1448 C 
(0.10), ~ 0.05. 

L 1527 is a special case. The Tnust (27 K) and especially 
the Lbol/F°3 ratio (5.0) suggest that it is a younger object, but 
its low N2D + /N2H + ratio and the fact that an embedded source 
is detected at wavelengths < 5/jm indicates that this source is 
more evolved (Froebrich 2005). Chemical models, which re- 
produce the observed abundance of long-chain unsaturated hy- 



drocarbons and cyanopolyynes, indicate a more evolved evolu- 
tionary stage of L 1527 as well (Hassel et a. 2008). From infra- 
red observations, Tobin et al. (2008) concluded that L 1527 
contains a massive protostellar disk, which is seen edge-on. 
Therefore L 1527 might be a more evolved object than indi- 
cated by continuum observations, due to the large obscuration 
of the central object. Hence, the N2D + /N2H + ratio seems to be 
even more reliable than continuum measurements to identify 
the youngest protostellar objects. 

In conclusion, the N2D + /N2H + ratio can be used to clearly 
indentify the youngest objects among the Class sources. For 
more evolved objects, the correlation of evolutionary stage and 
N2D + /N2H + ratio is not so clear, but an overall decrease of 
deuteration with evolutionary stage can be seen. 

5. Comparison with chemical- and radiative 
transfer models 

To explain the results discussed in the previous section, we ex- 
amined a model of a "typical" Class object rather than a spe- 
cific source. First, we describe the physical properties, e.g. the 
temperature and density profiles, of the model. Subsequently, 
we calculate the abundances of the important species and ap- 
ply a radiative transfer code to simulate the observations. 

5.1. Cloud structure 

The density structure of the model is that of a singular isother- 
mal sphere (Shu [T977T l: 

n(r) oc r~ 2 . 

Such density profiles have been found for Class protostars, 
too (Lonney et al. 2003). The inner boundary of the protostel- 
lar envelope is set to 100 AU. This limit is arbitrarily chosen, 
but at a distance of 220 pc, which is the distance of the Perseus 
molecular cloud, 100 AU corresponds to an angle of 0.46", 
which is much smaller than the angular resolution of our ob- 
servations. Therefore, the model results are not sensitive to the 
chosen inner boundary. The outer boundary is defined by the 
radius where the density is 10 4 cirr 3 . Kirk et al. (12006) esti- 
mated the external pressure on cores in the Perseus region and 
found values of logio(P e xt/k) between 5.5 and 6.0. Assuming a 
temperature between 15 and 30 K, this leads to a gas density 
between ~ 10 4 cnr 3 and ~ 5 ■ 10 4 cm" 3 . Around L 1448 IRS 3, 
Kwon et al. (|2006] > find a density of 6 - 9 • 10 3 cm -3 . Thus, 
10 4 crrT 3 can be considered an average value for the gas num- 
ber density surrounding Perseus cores. 

The total envelope mass is assumed to be of the order of one 
to a few solar masses. However, since the model results depend 
strongly on the density, the total mass varies for different mod- 
els (see below). For the temperature profile, we use a power 
law with an index of -0.4, which is the exponent for optically 
thin emission from protostellar envelopes, as found by Looney 
et al. (2003). They also showed that, with a density power-law 
index of 2, the radial temperature distribution can significantly 
deviate from the power-law at radii smaller than 100-200 AU. 
However, given that we used an inner boundary of 100 AU for 
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the envelope, such deviations do not noticeably influence our 
model. To avoid very low temperatures in the outskirts of the 
envelope, we assume a temperature profile like that found in 
prestellar cores (Crapsi et al. 120071) : 



T(r) = T out - 



Tout Ti n 

1 +(r/r ) L5 ' 



Here, we assume T out = 15 K, T; n =7 K, and ro=3633 AU. This 
radius corresponds to an angle of 26" at a distance of 220 pc, 
as found in L 1544. However, since we use this temperature 
profile only for the outer, low-density part of the envelope, the 
influence of these parameters on the model results is negligible. 
Density and temperature profiles are shown in Fig. [9] 
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Fig. 9. Density (dashed line) and temperature (solid lines) pro- 
files of the protostellar core models. The four temperature pro- 
files represent models that yield values for Tn U st of 20 K, 27 K, 
35 K and 45 K, respectively. 

Using this cloud structure, we calculate the continuum flux 
at 60 fim, 100 //m, 850 //m and 1.2 mm from the formula: 



(r). B v (T(r)) ■*,./! dr 



where n is the density, B(T) is the Planck function, and /j is the 
mass of the gas molecules. For the opacity per unit mass, k v , 
we use again 0.005 g/cm 2 at 1.2 mm, which assumes optically 
thin continuum emission. Furthermore, we assumed an emis- 
sivity index of f3 = 2, which is a typical value for molecular 
clouds (Ward-Thompson et al. 120021 ) and not the /3 — 1 shown 
by previous fits. We used /3 — 2, because the cold (~ 10 K) 
outskirts of the envelope still emit a significant amount of radi- 
ation at mm wavelengths, but much less at 60 /mi. The strong 
temperature gradient in protostellar cores causes the spectrum 
integrated over the entire core to be flattened. Thus, the emis- 
sivity index derived from a single temperature fit is underes- 
timated. The values for the continuum emission obtained by 
this method are then fitted, following the procedure described 
in section 14.21 and the absolute temperature scale is varied to 
obtain models with dust temperatures of 20 K, 27 K, 35 K and 
45 K, respectively, which covers the whole range of observed 
values. 



5.2. Chemical and radiative transfer model 

The chemical code used in this study is similar to the one orig- 
inally described by Caselli et al. 2002b but updated following 
Caselli et al. 2008 We assume a spherically symmetric cloud, 
and we follow the time-dependent freeze-out of CO and N2, as 
well as thermal and non-thermal desorption due to cosmic-ray 
impulsive heating, neglecting dynamical effects. All the deuter- 
ated forms of H^ are included in the model, and the abundances 
of the ionic species such as HCO + , N2H + and their deuterated 
couterparts are given by the steady state chemical equations 
using the instantaneous abundaces of the neutral species (see 
Caselli et al. 2002b for details). The rate coefficients have been 
adopted from the UMIST database for astrochemistry, avail- 
able at http://www.udfa.net For the proton-deuteron exchange 
reaction (such as H^ + HD — » H2D + + H2), we used both "stan- 
dard rates" (e.g. as used in Roberts et al. 2004 Ceccarelli & 
Dominik 2005) as well as the rates more recently measured 
by Gerlich et al. (2002), which appear to better fit our data 
(see below). No atomic oxygen is included in the code, given 
the large uncertainties in its abundance in cold and dense re- 
gions (e.g. Bergin & Snell l2002l see also discussion in section 
3.2.3 of Caselli et al. [2002b ]). The cosmic-ray ionization rate 
has been fixed at £=3xl0~ 17 s _1 , the average value measured 
toward high-mass star forming regions by van der Tak & van 
Dishoeck (2000). The CO and N2 binding energies are assumed 
to be 1 100 K and 982.3 K, respectively, following the prescrip- 
tions of Oberg et al. 2005 and Schnee et al. 2007 (see their sec- 
tion 6.2). The model is run for 10 6 yr, an average age for a dense 
core. After 10 years, the abundances of the various molecules 
have nearly reached chemical equilibrium, and thus the exact 
evolution time is not a critical parameter (see section IB31 for 
discussion). 

The results from the chemical model are used as input pa- 
rameters for a radiative transfer model, a non-LTE Monte-Carlo 
code developed by Pagani et al. (2007 ). This code, which is 
developed from the Bernes code (Bernes [19791 ), calculates the 
ID radiative transfer of both N2H + and N2D + and includes line 
overlap of the hyperfine transitions. Furthermore, recently pub- 
lished collisional coefficients between individual hyperfine lev- 
els have been implement (Daniel et al. 2005). For the convolu- 
tion with the beam, we assumed a Gaussian beam shape with 
HPBW as listed in Tab. [2] The distance to the cloud is set to 
220 pc, i.e. the distance of the Perseus cloud. From the spectra 
that we obtained from the Monte-Carlo-code, we determined 
the N2D + /N2H + ratio using the CTEX procedure, just as we did 
for the observed spectra (see section |4~T1 ). Hence, the compar- 
ison between the model and the observations is consistent, be- 
cause possible systematic errors, introduced by the assumption 
of a constant excitation temperature, are taken into account. 



5.3. Comparison between models and observations 

We compute the N2D + /N2H + ratios for four models with fit- 
ted dust temperatures of 20 K, 27 K, 35 K and 45 K, respec- 
tively. The deuterium fractionation in these models ranges from 
0.23 at 20 K to 0.03 at 45 K. These values fit the observations 
very well. We derive an H2 density at a radius of 100 AU of 
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2.77 ■ 10 8 crrT 3 , which corresponds to a total envelope mass of 
4.6 Mq. The modelled correlation between the deuterium frac- 
tionation and the CO depletion factor (f D (CO)) fits the obser- 
vations quite well, too, but in the cooler objects (Tpust ~ 20 K), 
the modelled CO freeze-out is very sensitive to the tempera- 
ture. Thus, a marginal temperature variation leads to a signifi- 
cant change of the CO depletion factor, and therefore the qual- 
ity of the fit depends strongly on the model used. In Fig. [10] a 
comparison of the models with the observations of the Perseus 
sources is shown. Results from the models are also listed in 
Table [TOl The best agreement with observations is found when 
using the reaction rates measured by Gerlich et al. 2002 for the 
proton-deuteron exchange reaction. Using the "standard rates", 
we could only fit the observed N2D + /N2H + ratios by increasing 
the total envelope mass to 25 M , a value that is clearly larger 
than expected for low mass protostars. 
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Fig. 10. Upper panel:N2D + /N2H + ratio versus dust tempera- 
ture. Lower panel: N2D + /N2H + ratio and the CO depletion 
factor Fd (CO). The open circles are the observed values for 
the sub-sample of cores in the Perseus cloud; the crosses are 
the model results. The N2D + /N2H + ratio of the models was de- 
termined using the CTEX method. 



To check that the assumption of a constant excitation tem- 
perature does not give a significantly different result from the 
T~ 04 profile, we calculated beam averaged column densities 
from the model abundance profiles as well. The N2D + /N2H + 
ratios determined from the beam averaged column densities are 
the same as before, to within 20%. The estimated error of the 
observed N2D + /N2H + ratio is also approximately 20%, so the 
fit is good, independently of the method used. The trend of a 
decreasing deuterium fraction with increasing To U st is true for 
both the ratio determined by the CTEX method and the ratio of 
the beam averaged column density. 



Table 10. Model results 



TdusI 


N 2 D + /N 2 H+ (CTEX) 


HC07DCO+ 


fo(CO) 


20 K 


0.23 


0.09 


4.73 


27 K 


0.14 


0.05 


2.55 


35 K 


0.07 


0.03 


1.97 


45 K 


0.03 


0.02 


1.67 



The N2D + J=3-2/J=2-l intensity ratios of the models range 
from 0.44 to 1.10, which is comparable to the observed ra- 
tios. For the J=2-l/J=l-0 ratios, the models give values a fac- 
tor of two higher than what is observed. However, these ratios 
are very sensitive to temperature, density and abundance varia- 
tions. Taking the simplicity of our model into account, we con- 
sider these values in agreement with the observations. 

Besides N 2 D + and N 2 H + , the abundances of HCO + and 
DCO + are also calculated by the chemical model. The 
DCO + /HCO + ratio is significantly lower than the N 2 D + /N 2 H + 
ratio and ranges from 0.02 to 0.08 in the models with dust 
temperatures of 45 K and 21 K, respectively. The HCO + and 
DCO + column densities are averaged over a 21" beam, and 
thus comparable to the observations conducted by J0rgensen et 
al. ([20041 . The values for the DCO + /HCO + ratio reported by 
them are between 0.004 and 0.048 and thus a little lower than 
the modelled ratios (see Tab.[T0l&l9ll. 

5.4. Chemical stratification within the core 

The temperature gradient, caused by the heating of the central 
object, leads to strong variations of the abundances of many 
molecules and molecular ions. The more immediate cause for 
these variations is the freeze-out of CO, which occurs at a 
temperature of ~ 20 K. The absence of CO in the gas phase 
has a big influence on the abundance of many other species, 
since many molecular ions, e.g. N2FU, are partially destroyed 
by reactions with CO. In addition, the deuterium fractionation 
is enhanced in regions where CO is depleted (see section [TJ. 
Therefore, the protostellar envelope can be roughly divided into 
three zones (Fig. fTTT i. In the inner zone, the temperature has 
increased above ~ 20 K, CO evaporates from the grains and, 
consequently, the abundance of N2FF as well as the deuterium 
fractionation are relatively low in this region. As soon as the 
temperature drops below the critical value of ~ 20 K, the CO 
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abundance drops drastically. At the same time, the N2H + abun- 
dance rises by 1-2 orders of magnitude. The N2D + /N2H + ratio 
rises as well, but its maximum is at a larger radius than the 
CO abundance minimum. This is because the temperature at 
the minimum of the CO abundance is ~ 20 K, where the de- 
struction of H 2 D + via the reaction H 2 D + + H 2 -> + HD 
is not neglible. At even larger radii, the degree of CO deple- 
tion decreases as the gas density decreases, which causes the 
N 2 H + abundance to stay more or less constant, and the deu- 
terium fractionation to drop again. 



The reason for the different ratios of the deuterated and 
non-deuterated forms of HCO + and N 2 H + is that the deuterium 
fractionation is strongest in places where CO, a progenitor of 
HCO + & DCO + , is mostly depleted. Thus, contrary to N 2 H + , 
HCO + is also relatively strong in the warm and dense innner 
part of the envelope, and the ratio of the column densities of 
DCO + and HCO + is significantly lower than the N 2 D + /N 2 H + 
ratio. This can be seen in Fig. [T2] in which the abundance of 
the molecules and molecular ions are shown as a function of 
radius. 




r [AU] 



Fig. 11. The relative abundance of CO (upper panel), N 2 H + 
(middle panel) and the N 2 D + /N 2 H + ratio are shown as a func- 
tion of radius. It can be clearly seen that the radius at which CO 
starts to deplete, as well as the peak of the deuterium fraction- 
ation, shifts outwards with increasing temperature. 
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Fig. 12. The abundances of CO, H+, H 2 D + , HCO + , DCO + , 
N 2 H + , and N 2 D + relative to H 2 are shown as a function of ra- 
dius. The thick lines show the abundance profile of the proto- 
nated species. The thin lines show their singly-deuterated coun- 
terparts. These profiles apply to the model with To us t = 27 K. 

The radial abundance profiles we found with our static 
model are in good agreement with the results of Lee et 
al. d20051 l. who calculated a self-consistent dynamical and 
chemical model of the entire star forming process. In their 
model, the deuterated species are not calculated, but the abun- 
dance profiles of CO, N 2 H + and HCO + at a time ~ 10 5 yr af- 
ter the collapse starts are similar to the profiles we found (see 
their Fig. 11). This similarity indicates that the assumption of a 
chemical equilibrium is valid at least for the simple species. 
In a more recent paper, Aikawa et al. (|2008) get somewhat 
shorter timescales for the protostellar evolution: 9.3 • 10 4 years 
after the beginning of the collapse, the CO desorption radius 
is ~ 2000 AU, which is comparable to the results of our 45 K 
model. But even for these shorter time scales, the assumption 
of chemical equilibrium is still valid (see next section). 



During the Class phase, the temperature of the central 
object increases, and so does the temperature of the enve- 
lope. Consequently, the radius at which CO freezes out onto 
dust grains shifts further outwards, and the fraction of the 
cloud in which the abundance of deuterated molecules is en- 
hanced gets smaller with time. Furthermore, the mean density 
of the deuterium fractionated gas decreases, thus the mean local 
N 2 D + /N 2 H + ratio declines as well. This leads to the decrease 
of the N 2 D + /N 2 H + column density ratio, which we observe in 
our sample. 



5.5. Sensitivity of the Model 

The evolution times of the models are the least constrained 
values. The estimated ages of Class protostars range from 
~ 10 4 to a few times 10 5 years after the collapse has started 
(Froebrich 2005). The ages vary by a factor of 10 depending 
which model is used. However, the crucial time for the deu- 
terium chemistry is the time since the CO started to freeze out, 
because the freeze-out of CO is the slowest reaction in the net- 
work. Thus, the entire prestellar phase has to be taken into ac- 
count. After a time of 10 6 years, the abundances have already 
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reached equilibrium. Reducing the time to 10 5 years changes 
the N2D + /N2H + ratio only marginally. The difference between 
the ratios after 10 6 years and 10 5 years is less than 0.01. In 
models with a evolution time of 10 4 years, the abundances 
are no longer in equilibrium. In the coldest protostar (20 K), 
the deuterium fractionation is already high (0.26), but for the 
warmer objects, the N2D + /N2H + ratio is a factor ~ 2 lower than 
in chemical equilibrium. Although already quite low, these val- 
ues still agree with the observations. Evolution times longer 
than 10 7 years do not change the results, because equilibrium 
is already reached. 

Another important parameter is the cosmic-ray ionisation 
rate. In our model, we used f = 3 ■ 10 17 s _1 (van der Tak & 
van Dishoeck l2000l i. A f of 6 • 10 17 s _I increases the N 2 H + 
abundance by a factor of 1.5, but the deuterium fractionation is 
lowered by 40%. A cosmic -ray ionization rate of 1.5 ■ 10 17 s" 1 
increases the N 2 D + /N 2 H + ratio by 30%. 

The central density of the core, and consequently the enve- 
lope mass, has a similar effect on the deuterium fractionation 
as In our models, the density is varied to match the observed 
N2D + /N2H + ratio best. The only constraint we put on the cen- 
tral density is that the resulting envelope mass should be of the 
order of a few solar masses. For our best fit model, we derive 
an envelope mass of 4.6 M . Decreasing the density by a factor 
of a half, the deuterium fractionation decreases by 40%. With 
a density twice as high as our best fit model, the N2D + /N2H + 
increases by a factor of 1 .2. 

The last crucial parameter we investigated is the binding 
energy of CO and N2 on dust grains. In our model, we used 
binding energies of 1100 K and 982.3 K for CO and N 2 , re- 
spectively, which are the binding energies on mixed CO-H2O 
ice. We also calculated the models using the values for Eb on 
pure CO ice (855 and 790 for CO and N2, respectively, Oberg 
et al. l2005b . With these values, the resulting N 2 D + /N 2 H + ratios 
drop to < 0.04 for all four temperatures. 

6. Summary 

We observed 20 Class protostars in N 2 H + 1-0, N 2 D + 1-0, 
2-1 and 3-2, C ls O 1-0 and HCO + 3-2. The integrated inten- 
sities of the N2H + 1-0 lines in our sample are comparable to 
those found in prestellar cores (Crapsi et al. [2005 ). The op- 
tical depths are typically 67% lower in the Class sources, 
but because the excitation temperature is 2-3 K higher and the 
line width is about 2.5 time lager, the N2H + column densities 
(~ 10 13 crrT 2 ) are also comparable. The presence of a protostar 
also affects the kinematics of the cores. N2H + lines are sig- 
nificantly wider in Class sources than in prestellar cores (on 
average, 0.61 km/s and 0.26 km/s, respectively). A comparison 
of our HCO + observations with observations conducted with a 
-three times larger beam (Gregersen et al. 2000) leads us to the 
conclusion that HCO + 3-2 stems mainly from a compact region 
of about 4000 AU in size. 

The dependence of the N2D + /N2H + ratio on dust temper- 
ature, CO depletion factor, and the Lbol/F^ 3 pc ratio clearly 
indicates a close relation between deuterium fractionation and 
evolutionary stage of a Class protostar. For the sub-sample 
of sources in the Perseus cloud, the correlation is striking. This 



might be an indication for the influence of the environment on 
the deuterium fractionation. 

The correlation of the deuterium fractionation and the 
CO depletion factor looks qualitatively like the one found 
in prestellar cores (Crapsi et al. 120051 1. However, whereas in 
prestellar cores, the N2D + /N2H + ratio stays low until fo(CO) ~ 
10, in protostars, the ratio starts to rise at fo(CO) ~ 3. This dif- 
ference occurs because CO is highly depleted at the centre of 
prestellar cores, but not at the center of protostellar cores. 

There is a weak correlation of N2D + line width with 
N2D + /N2H + ratio, but line broadening by multiple sources 
makes this conclusion suspect. A clear correlation with other 
kinematical tracers, especially with 6v, is found. Sources with a 

high N2D + /N2H + ratio clearly show infall motion (6v 0.4). 

The lower the deuterium fractionation gets, the more 5v in- 
creases. This might reflect the fact that outflow activity in- 
creases with evolutionary stage (Richer et al. 2000). 

A model with a power law density and temperature profile 
reproduces the correlation between the N2D + /N2H + ratio and 
both Toust and the CO depletion factor very well. In these mod- 
els, the protostellar envelopes are chemically stratified. The in- 
ner part is too warm for CO to freeze out, so that the N2FF 
abundance and the deuterium fractionation are low. At radii 
where the temperature drops below ~ 20 K, CO is mostly de- 
pleted and the N2D + /N2H + ratio increases by several orders of 
magnitude. At even larger radii, the lower density decreases 
the degree of CO depletion, and thus the N2D + /N2H + ratio also 
decreases. 

The results of this work, in combination with the results 
of Crapsi et al. d2005t show that the deuterium fractionation 
is highest at the moment of collapse. The primary use of the 
deuterium fractionation of N2H + is to identify very young pro- 
tostellar objects, because the decline of the N2DVN2I-F ratio 
is rather steep. 

Comparisons of N2FF fractionation with the deuterium 
fractionation of HCO + and NH3 show clear differences, 
which are due either to the influence of freeze-out on dust 
grains or differences in chemical reactions or both. The low 
DCO + /HCO + ratio observed by J0rgensen et al. (120041 1 is 
a prediction of our model. In young, i.e. cool objects, the 
NH2D/NH3 ratio, observed by Hatchell (120031) is comparable 
to the deuterium fractionation of N2H + . As the temperature 
increases, the N2D + /N2H + ratio falls off much faster than the 
NH2D/NH3 ratio, which stays above 0.15 in all observed ob- 
jects. This cloud be due to desorption of deuterium enhanced 
ammonia from the dust grains. 
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